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I . INTRODUCTION 

The Penn State Finite Difference Time Domain Electromagnetic 
Scattering Code Version C is a three dimensional numerical 
electromagnetic scattering code based upon the Finite Difference 
Time Domain Technique (FDTD) . The supplied version of the code 
is one version of our current three dimensional FDTD code set. 
This manual provides a description of the code and corresponding 
results for several scattering problems. The manual is organized 
into fourteen sections: introduction, description of the FDTD 

method, operation, resource requirements. Version C code 
capabilities, a brief description of the default scattering 
geometry, a brief description of each subroutine, a description 
of the include file (COMMONC. FOR) , a section briefly discussing 
Radar Cross Section (RCS) computations, a section discussing some 
scattering results, a sample problem setup section, a new problem 
checklist, references and figure titles. 

II. FDTD METHOD 

The Finite Difference Time Domain (FDTD) technique models 
transient electromagnetic scattering and interactions with 
objects of arbitrary shape and/or material composition. The 
technique was first proposed by Yee [1] for isotropic, non- 
dispersive materials in 1966; and has matured within the past 
twenty years into a robust and efficient computational method. 

The present FDTD technique is capabable of transient 
electromagnetic interactions with objects of arbitrary and 
complicated geometrical shape and material composition over a 
large band of frequencies. This technique has recently been 
extended to include dispersive dielectric materials, chiral 
materials and plasmas. 

In the FDTD method. Maxwell's curl equations are discretized 
in time and space and all derivatives (temporal and spatial) are 
approximated by central differences. The electric and magnetic 
fields are interleaved in space and time and are updated in a 
second-order accurate leapfrog scheme. The computational space 
is divided into cells with the electric fields located on the 
edges and the magnetic fields on the faces (see Figure 1) . FDTD 
objects are defined by specifying dielectric and/or magnetic 
material parameters at electric and/or magnetic field locations. 

Two basic implementations of the FDTD method are widely used 
for electromagnetic analysis: total field formalism and scattered 
field formalism. In the total field formalism, the electric and 
magnetic field are updated based upon the material type present 
at each spatial location. In the scattered field formalism, the 
incident waveform is defined analytically and the scattered field 
is coupled to the incident field through the different material 
types. For the incident field, any waveform, angle of incidence 
and polarization is possible. The separation of the incident and 



scattered fields conveniently allows an absorbing boundary to be 
employed at the extremities of the discretized problem space to 
absorb the scattered fields. 
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This code is a scattered field code, and the total E and H 
fields may be found by combining the incident and scattered 
fields. Any type of field quantity (incident, scattered, or 
total) , Poynting vector or current are available anywhere within 
the computational space. These fields, incident, scattered and 
total, may be found within, on or about the interaction object 
placed in the problem space. By using a near to far field 
transformation, far fields can be determined from the near fields 
within the problem space thereby affording radiation patterns and 
RCS values. The accuracy of these calculations is typically 
within a dB of analytic solutions for dielectric and magnetic 
sphere scattering. Further improvements are expected as better 
absorbing boundary conditions are developed and incorporated. 

III. OPERATION 

Typically, a truncated Gaussian incident waveform is used to 
excite the system being modeled, however certain code versions 
also provide a smooth cosine waveform for convenience in modeling 
dispersive materials. The interaction object is defined in the 
discretized problem space with arrays at each cell location 
created by the discretization. All three dielectric material 
types for E field components within a cell can be individually 
specified by the arrays IDONE(I, J,K) , IDTWO(I, J,K) , 

IDTHRE(I, J,K) . This models arbitrary dielectric materials with ju 
=./V an obvious extension to six arrays, magnetic materials 

with n f Hq can be modeled. 

Scattering occurs when the incident wave, marched forward in 
time in small steps set by the Courant stability condition, 
reaches the interaction object. Here a scattered wave must 
appear along with the incident wave so that the Maxwell equations 
are satisfied. If the material is a perfectly conductive metal 
then only the well known boundary condition 

scat inc . „ . 

®tan = “E tan (1) 

must be satisfied. For a nondispersive dielectric the 
requirement is that the total field must satisfy the Maxwell 
equations in the material: 


Vx E = Vx ( E lnc + E scat ) =- — 


1 8H* 


1 <9 (H + H scat ) 


M 0 dt 
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VxH tot =Vx (H ,nc +H scat ) = c 


dE 


tot 


at 


+ CT E 


tot 


(3) 


= e 


a ( E inc + E scat ) 


at 


+ a (E + E ) 


(4) 


Additionally the incident wave, defined as moving unimpeded 
through a vacuum in the problem space, satisfies everywhere in 
the problem the Maxwell equations for free space 


Vx E = - — 


i aH’ 


dt 


(5) 


VxH ,nc = e 


aE' 


at 


( 6 ) 


Subtracting the second set of equations from the first yields the 
Maxwell equations governing the scattered fields in the material: 


nScat i aH scat 


Vx E sca — 


< 3 t 


(7) 


VxH scat = (e -e ) 


aE 


inc 


at 


+ aE ,nc + e 


dE 


,scat 


at 


+ a E 


scat 


( 8 ) 


Outside the material this simplifies to: 

1 aH scat 


Vx E 


scat 


M 0 dt 


(9) 


VxH scat =c 


aE scat 

at 


( 10 ) 


Magnetic materials, dispersive effects, non-linearities, 
etc., are further generalizations of the above approach. Based 
on the value of the material type, the subroutines for 
calculating scattered E and H field components branch to the 
appropriate expression for that scattered field component and 
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that component is advanced in time according to the selected 
algorithm. As many materials can be modeled as desired, the 
number eguals the dimension selected for the flags. If materials 
with behavior different from those described above must be 
modeled, then after the appropriate algorithm is found, the 
code's branching structure allows easy incorporation of the new 
behavior. 

IV. RESOURCE REQUIREMENTS 

The number of cells the problem space is divided into times 
the six components per cell set the problem space storage 
requirements 

Storage = NC x 6 components/cell x 4 bytes /component (11) 
and the computational cost 

Operations = NC x 6 comp/cell x 10 ops/component x N (12) 
where N is the number of time steps desired. 

N typically is on the order of ten times the number of cells 
on one side of the problem space. More precisely for cubical 

cells it takes ^ time steps to traverse a single cell when the 
time step is set by the Courant stability condition 

Ax 

At = Ax = cell size dimension (13) 

/3c 


The condition on N is then that 

N - 10X (/rj, Nc'= - number ° ellS ° n a Side < 14 > 

of the problem space 

The earliest aircraft modeling using FDTD with approximately 30 
cells on a side required approximately 500 time steps. For more 
recent modeling with approximately 100 cells on a side, 2000 or 
more time steps are used. 

For (100 cell) 3 problem spaces, 24 MBytes of memory are 
required to store the fields. Problems on the order of this size 
have been run on a Silicon Graphics 4D 220 with 32 MBytes of 
memory, IBM RISC 6000, an Intel 486 based machine, and VAX 
11/785. Storage is only a problem as in the case of the 486 
where only 16 MBytes of memory was available. This limited the 
problem space size to approximately (80 cells) 3 . 
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For (100 cell) problems with approximately 2000 time steps, 
there is a total of 120 x 10 operations to perform. The speeds 
of the previously mentioned machines are 24 MFLOPs (4 processor 
upgraded version), 10 MFLOPS, 1.5 MFLOPS, and 0.2 MFLOPs. The 
run times are then 5^ x 10 seconds, 12 x 1CT seconds, 80 x 1CT 
seconds and 600 x 10 seconds, respectively. In hours the times 
are 1.4, 3.3, 22.2 and 167 hours. Problems of this size are 
possible on all but the last machine and can in fact be performed 
on a personal computer (486) if one day turnarounds are 
permissible. 

V. VERSION C CODE CAPABILITIES 

The Penn State University FDTD Electromagnetic Scattering 
Code Version C has the following capabilities: 

1) Ability to model lossy dielectric and perfectly conducting 
scatterers. 

2) Ability to model lossy magnetic scatterers. 

3) First and second order outer radiation boundary condition 
(ORBC) operating on the electric fields for dielectric 
scatterers . 

4) First and second order ORBC operating on the magnetic fields 
for magnetic scatterers. 

5) Near to far zone transformation capability to obtain far zone 
scattered fields. 

6) Gaussian and smooth cosine incident waveforms with arbitrary 
incidence angles. 

7) Near zone field, current or power sampling capability. 

8) Companion code for computing Radar Cross Section (RCS) . 

VI. DEFAULT SCATTERING GEOMETRY 

The code as delivered is set up to calculate the far zone 
backscatter fields for a 20 4 cm radius, lossy magnetic sphere with 
parameters e=e 0 , n=4 n 0 and a=2.8413E+4. The problem space size 
is 66 by 66 by 66 cells in the x, y and z directions, the cells 
0»01 cm cubes, and the incident waveform is a 0— polarized 
Gaussian pulse with incidence angles of 0=22.5 and 0=22.5 
^®9 rees * The output data files are included as a reference along 
with a code (RCS3D . FOR) for computing the frequency domain RCS 
using these output data files. The ORBC is the second order 
absorbing boundary condition set forth by Mur (2]. 
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VII. SUBROUTINE DESCRIPTION 

In the description for each subroutine, an asterisk (*) will 
be placed by the subroutine name if that particular subroutine is 
normally modified when defining a scattering problem. 

MAIN ROUTINE 

The main routine in the program contains the calls for all 
necessary subroutines to initialize the problem space and 
scattering object (s) and for the incident waveform, far zone 
transformation, field update subroutines, outer radiation 
boundary conditions and field sampling. 

The main routine begins with the include statement and then 
appropriate data files are opened, and subroutines ZERO, BUILD 
and SETUP are called to initialize variables and/or arrays, build 
the object (s) and initialize the incident waveform and 
miscellaneous parameters, respectively. Subroutine SETFZ is 
to intialize parameters for the near to far zone 
"transformation if far zone fields are desired. 

The main loop is entered next, where all of the primary 
field computations and data saving takes place. During each time 
step cycle, the EXSFLD, EYSFLD, and EZSFLD subroutines are called 
to update the x, y, and z components of the scattered electric 
field. The six electric field outer radiation boundary 
conditions (RADE??) are called next to absorb any outgoing 
scattered fields for perfectly conducting or lossy dielectric 
scatterers. Time is then advanced 1/2 time step according to the 
Yee algorithm and then the HXSFLD, HYSFLD, AND HZSFLD subroutines 
are called to update the x, y, and z components of scattered 
magnetic field. The six magnetic field outer radiation boundary 
conditions (RADH??) are called next to absorb any outgoing 
scattered fields for lossy magnetic scatterers. Time is then 
advanced another 1/2 step and then either near zone fields are 
sampled and written to disk in DATSAV, and/or the near zone to 
far zone vector potentials are updated in SAVFZ. The parameter 
NZFZ (described later) in the common file defines the type of 
output fields desired. 

After execution of all time steps in the main field update 
loop, subroutine FAROUT is called if far zone fields are desired 
to compute the far zone fields and write them to disk. At this 
point, the execution is complete. 

SUBROUTINE SETFZ 

This subroutine initializes the necessary parameters 
required for far zone field computations. The code as furnished 
computes backscatter far zone fields and can compute bistatic far 
zone fields for one scattering angle (i.e. one 0 and 0 angle) . 
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Refer to reference [3] for a complete description of the near to 
far zone transformation. Other versions of this subroutine 
provide for multiple bistatic angles. 

SUBROUTINE SAVFZ 


This subroutine updates the near zone to far zone vector 
potentials. 

SUBROUTINE FAROUT 

This subroutine changes the near zone to far zone vector 
potentials to far zone electric field 0 and <p components and 
writes them to disk. 

SUBROUTINE BUILD * 


This subroutine "builds" the scattering object (s) by 
initializing the I DONE, IDTWO, IDTHRE, IDFOR, IDFIV and IDSIX 
arrays. The I DONE- IDTHRE arrays are for specifying perfectly 
conducting and lossy dielectric materials. The IDFOR-IDSIX 
arrays are for lossy magnetic materials. The reason for the 
separate arrays is so the user can independently control the 
exact placement of dielectric and magnetic material in the Yee 
cells. Refer to Figure l for a diagram of the basic Yee cell. 

For example, setting an element of the I DONE array at some I,J,K 
location is actually locating dielectric material at a cell edge 
center location is I+0.5,J,K. Setting an element of the 
IDFOR array at some I,J,K location is actually locating magnetic 
material perpendicular to a cell face whose center location is 
I, J+0. 5 , K+0 . 5 , or equivalently, locating magnetic material at an 
edge on the dual H field mesh. The difference between the IDONE 
and IDFOR array locations is a direct result of the field offsets 
m the Yee cell (see Figure 1). Thus, materials with diagonal 
permitti^y anc */° r diagonal permeability tensors can be modeled. 
The default material type for all ID??? arrays is 0, or free 
space. By initializing these arrays to values other than 0, the 
user is defining an object by determining what material types are 
present at each spatial location. Other material types available 
for I DONE- IDTHRE are 1 for perfectly conducting objects, and 2-9 
for lossy dielectrics. I DONE- IDTHRE are normally set to 0 for 
magnetic scatterers. Other material types available for IDFOR- 
IDSIX are 10-19 for lossy magnetic materials. IDFOR-IDSIX are 
normally set to 0 for perfectly conducting or dielectric 
scatterers. If the user wants a material with both dielectric 
and magnetic properties (i.e. permittivity other than e n for 
magnetic materials, and permeability other than n Q for dielectric 
materials), then he/ she must define IDONE-IDSIX for that 
particular material. This subroutine also has a section that 
checks the ID??? arrays to determine if legal material types have 
been defined throughout thf problem space. The actual material 
parameters (e, /x, a, and a ) are defined in subroutine SETUP. 



The default geometry is a 20 cm radius, lossy magnetic sphere. 

The user must be careful that his/her object created in the 
BUILD subroutine is properly formed. There is not a direct 
one-to-one correspondence between the dielectric and magnetic 
ID??? arrays. However, one can define a correspondence, so that 
code used to generate a dielectric object can be modified to 
generate a magnetic object. 

To see this consider that we have set the permittivity at 
cell locations corresponding to 

EX ( I , J , K) , EY ( I , J , K) , EZ ( I , J , K) 

using the IDONE, IDTWO, and IDTHRE arrays respectively. This 
determines one corner of a dielectric cube. If we wish to define 
the corner of a corresponding magnetic cube, offset 1/2 cell in 
the x, y, z directions, we would set the locations of the 
magnetic fields 


HX(I+1,J,K), HY(I,J+1,K), HZ ( I , J, K+l) 

as magnetic material using the IDFOR, IDFIV, and IDSIX arrays. 

This example indicates the following correspondence between 
BUILDing dielectric and magnetic objects: 


Dielectric Magnetic 

Object Object 


IDONE ( I , J, K) = IDFOR ( 1+1 , J, K) 

IDTWO ( I , J , K) = IDFIV ( I , J+l , K) 

IDTHRE ( I , J , K) = IDSIX ( I , J , K+l) 

T° illustrate this for a somewhat more complicated case, consider 
an example of a 3x3 cell dielectric plate located in the XY plane 
at K=5. The plate can be generated with the following FORTRAN 
lines: 


7 


DO 20 J=4 , 7 6 

DO 10 1=4,7 

IF ( I . NE . 7 ) IDONE ( I , J, K) =2 J 

IF ( J. NE . 7 ) IDTWO ( I , J, K) =2 * 5 

10 CONTINUE | 

20 CONTINUE j 

! 4 

j 4 5 6 7 


PLATE 


f 
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If the same FORTRAN lines were used to try to generate a magnetic 
plate, the object generated would actually be unconnected: 


7 


PLATE 


DO 20 J=4 , 7 
DO 10 1=4,7 

IF ( I . NE . 7 ) IDFOR (I , J, K) =11 
IF ( J . NE . 7 ) IDFIV ( I , J , K) =11 
10 CONTINUE 
20 CONTINUE 



The correct way to build the magnetic plate is with the followinq 
FORTRAN code : 


DO 20 J=4 , 7 
DO 10 1=4,7 

IF ( I . NE . 4 ) IDFOR (I,J,K)=11 
IF ( J . NE . 4 ) IDFIV (I,J,K)=11 
10 CONTINUE 
20 CONTINUE 

This is equivalent to the correspondence relationship given 
above. See comments in the BUILD subroutine for further 
explanation of the ID??? arrays. 

When it is important to place the object in the center of 
b problem space (to have lowest possible cross— pol scattering 
for symmetric objects), NX etc. should be odd for dielectrics and 
even for magnetics. This is due to the field locations in the 
Yee cell and also the placement of the E and H field absorbing 
boundary condition surfaces. 

If the object being modeled has curved surfaces, edges, etc. 
hat are at an angle to one or more of the coordinate axes, then 
hat shape must be approximately modeled by lines and faces in a 
stair-stepped" (or stair-cased) fashion. This stair-cased 
approximation introduces errors into computations at higher 
frequencies. Intuitively, the error becomes smaller as more 
cells are used to stair-case a particular object. Thus, the 
default Nickel Ferrite sphere scattering geometry is a stair- 
cased sphere. 

When the user's test object is dielectric it is apparently 
better to use the Mur E field boundary formulation. For a 
magnetic scattering object it seems better to use the Mur H field 
boundary formulation. The reason for this has to do with the 
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reactive part of the energy (near zone fields) reaching the 
absorbing boundary surface. These near zone fields differ for 
magnetic and dielectric objects, even for dual objects where the 
far fields are essentially identical. 

One final comment on building dual dielectric vs magnetic 
objects can be explained by considering the duality between 
magnetic and dielectric materials. In testing the code, it was 
helpful to test dual dielectric/magnetic objects since they 
scatter identically in the far zone. In specifying the material 
of a magnetic scatterer to be the dual of a dielectric, the 
following transformation must be applied: 

DIELECTRIC < — DUAL — > MAGNETIC 


E field 
H field 
€ 

6 (k) 

a 


< > H field 

< > -e field 

< > n 

< > e 

< > B (k) 


> <?• (M 0 /e 0 )=CT* 


The somewhat surprising entry in this table is the relationship 
between the dielectric and magnetic conductivities. The reason 
why the magnetic conductivity, a , for magnetic materials has to 
be scaled as indicated is that for duality to be applied the free 
space impedance must be inverted. Although this is not generally 
possible for actual problems, identical far zone scattering for 
dielectric apd magnetic scatterers can be achieved by this 
scaling of a as above and realizing that E and H scattered field 
incident . field ratios will not invert. This scaling of a is 
why conductivity for magnetic materials is usually not a 
dominating feature, and in fact is often neglected. 


SUBROUTINE DCUBE 


This subroutine builds cubes of dielectric material by 
defining four each of I DONE , IDTWO and IDTHRE components 
corresponding to one spatial cube of dielectric material. It can 
also be used to define thin (i.e. up to one cell thick) 
dielectric or perfectly conducting plates. Refer to comments 
within DCUBE for a description of the arguments and usage of the 
subroutine. This subroutine could be modified to build cubes 
and/or plates of magnetic materials by using a triple do-loop in 
BUILD (after calls to DCUBE) over coordinate indices I, J, K and 
applying the correspondence between the IDONE-IDTHRE and IDFOR- 
IDSIX arrays that was previously discussed. 

SUBROUTINE SETUP * 

This subroutine initializes many of the constants required 
for incident field definition, field update equations, outer 
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radiation boundary conditions and t material parameters. The 
material parameters e, fi, a and a are defined for each material 
type using the material arrays EPS, XMU, SIGMA and SIGMAC 
respectively. The array EPS is used for the total permittivity, 
XMU is used for the total permeability, SIGMA is used for the 
electric conductivity and SIGMAC is used for the magnetic 
conductivity (useful for running dual problems) . These arrays 
are initialized in SETUP to free space material parameters for 
all material types and then the user is required to modify these 
arrays for his/her scattering materials. Thus, for the lossy 
dielectric material type 2, the user must define EPS (2) and 
SIGMA (2) . The remainder of the subroutine computes constants 
used in field update equations and boundary conditions and writes 
the diagnostics file. 

SUBROUTINE EXSFLD 


fhi s subroutine updates all x components of scattered 
electric field at each time step except those on the outer 
boundaries of the problem space. IF statements based upon the 
I DONE array are used to determine the type of material present 
and the corresponding update equation to be used. These 

f field equations are based on the development given in 

[ 4 ] • 


SUBROUTINE EYSFLD 


This subroutine updates all y components of scattered 
electric field at each time step except those on the outer 
boundaries of the problem space. IF statements based upon the 
IDTWO array are used to determine the type of material present 
and the corresponding update equation to be used. 

SUBROUTINE EZSFLD 


This subroutine updates all z components of scattered 
electric field at each time step except those on the outer 
boundaries of the problem space. IF statements based upon the 
IDTHRE array are used to determine the type of material present 
and the corresponding update equation to be used. 

SUBROUTINES RADEYX, RADEZX, RADEZY, RADEXY, RADEXZ and RADEYZ 

, These subroutines apply the outer radiation boundary 
conditions to the scattered electric field on the outer 
boundaries of the problem space for non— magnetic scatterers . The 
user controls selection of these routines through the parameter 
MAGNET (described later) . 
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SUBROUTINE HXSFLD 

This subroutine updates all x components of scattered 
magnetic field at each time step. IF statements based upon the 

IDFOR array are used to determine the type of material present 
and the corresponding update equation to be used. 

SUBROUTINE HYSFLD 

This subroutine updates all y components of scattered 
magnetic field at each time step. IF statements based upon the 
IDFIV array are used to determine the type of material present 
and the corresponding update equation to be used. 

SUBROUTINE HZSFLD 

This subroutine updates all z components of scattered 
magnetic field at each time step. IF statements based upon the 
IDSIX array are used to determine the type of material present 
and the corresponding update equation to be used. 

SUBROUTINES RADHXZ , RADHYX, RADHZY, RADHXY, RADHYZ and RADHZX 

These subroutines apply the outer radiation boundary 
conditions to the scattered magnetic field on the outer 
boundaries of the problem space for magnetic scatterers. The 
user controls selection of these routines through the parameter 
MAGNET (described later) . 

SUBROUTINE DATSAV * 

This subroutine samples near zone scattered field quantities 
and saves them to disk. This subroutine is where the quantities 
to be sampled and their spatial locations are to be specified and 
is only called if near zone fields only are desired or if both 
near and far zone fields are desired. Total field quantities can 
also be sampled. See comments within the subroutine for 
specifying sampled scattered and/or total field quantities. When 
sampling magnetic fields, remember the <St/2 time difference 
between E and H when writing the fields to disk. Sections of 
code within this subroutine determine if the sampled quantities 
and the spatial locations have been properly defined. 

FUNCTIONS EXI, EYI, EZI, HXI , HYI and HZI 

These functions are called to compute the x, y and z 
components of incident electric and magnetic field, respectively. 
The functional form of the incident field is contained in a 
separate function SOURCE. 
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FUNCTION SOURCE * 

. This function contains the functional form of the incident 
field. The code as furnished uses the Gaussian form of the 
incident field. An incident smooth cosine pulse is also 
available by uncommenting the required lines and commenting out 
the Gaussian pulse. Thus, this function need only be modified if 
the user changes the incident pulse from Gaussian to smooth 
cosine. A slight improvement in computing speed and 
vectorization may be achieved by moving this function inside each 
of the incident field functions EXI, EYI and so on. 

FUNCTIONS DEXI, DEYI , DEZI, DHXI , DHYI and DHZI 

These functions are called to compute the x, y and z 
components of the time derivative of incident electric and 
magnetic field, respectively. The functional form of the 
incident field is contained in a separate function DSRCE. 

FUNCTION DSRCE * 


This function contains the functional form of the time 
derivative of the incident field. The code as furnished uses the 
time derivative of the Gaussian form of the incident field. A 
smooth cosine pulse time derivative is also available by 
uncommenting the required lines and commenting out the Gaussian 
Thus, the function need only be modified if the user 
changes from the Gaussian to smooth cosine pulse. Again, a 
slight improvement in computing speed and vectorization may be 
achieved by moving this function inside each of the time 
derivative incident field functions DEXI, DEYI and so on. 

SUBROUTINE ZERO 

This subroutine initializes various arrays and variables to 

zero. 


VIII. INCLUDE FILE DESCRIPTION (COMMONC. FOR) * 

The include file, COMMONC. FOR, contains all of the arrays 
and variables that are shared among the different subroutines. 
This file will require the most modifications when defining 
scattering problems. A description of the parameters that are 
normally modified follows. 

The parameters NX, NY and NZ specify the size of the problem 
space in cells in the x, y and z directions respectively. For 
problems where it is crucial to center the object within the 
problem space, then NX, NY and NZ should be odd for dielectric 
scatterers and even for magnetic scatterers. The parameter NTEST 
defines; the number of near zone quantities to be sampled and NZFZ 
defines the field output format. Set NZFZ=0 for near zone fields 
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only, NZFZ=1 for far zone fields only and NZFZ=2 for both near 
and far zone fields. Parameter MAGNET is used to define magnetic 
scatterers and it controls the choice of RADE?? versus RADH?? 
absorbing boundary subroutines. It is set to 0 for dielectric 
scatterers and to 1 for magnetic scatterers. Parameter NSTOP 
defines the maximum number of time steps. DELX, DELY, and DELZ 
(in meters) define the cell size in the x, y and z directions 
respectively. The 0 and 0 incidence angles (in degrees) are 
defined by THINC and PHINC respectively and the polarization is 
defined by ETHINC and EPHINC. ETHINC=1.0, EPHINC=0.0 for 0- 
polarized incident field and ETHINC=0.0, EPHINC=1.0 for 0- 
polarized incident fields. Parameters AMP and BETA define the 
maximum amplitude and the e temporal width of the incident 
pulse respectively. BETA automatically adjusts when the cell 
size is changed and normally should not be changed by the user. 
The far zone scattering angles are defined by THETFZ and PHIFZ. 
The code as furnished performs backscatter computations, but 
these parameters could be modified for a bistatic computation. 

IX. RCS COMPUTATIONS 

A companion code, RCS3D.F0R, has been included to compute 
RCS versus frequency. It uses the file name of the FDTD far zone 
data (FZ0UT3D.DAT) and writes data files of far zone 
electric fields versus time (FZTIME.DAT) and RCS versus frequency 
(3DRCS.DAT). The RCS computations are performed up to the 10 
cell/A 0 frequency limit. Refer to comments within this code for 
further details. 

X. RESULTS 

As previously mentioned, the code as furnished models a 20 
cm radius, lossy magnetic sphere and computes backscatter far 
zone scattered fields at angles of 0=22.5 and 0=22.5 degrees. 
Results are included for the dual dielectric sphere and the 
default magnetic sphere. The material parameters for the 
dielectric dual of the magnetic material are e=4e 0 , cr=0.2 and 
M = M o* For these materials there are 5 cells per wavelength at 
approximately 3 . 0 GHz . 

Figures 2-3 show the co-polarized far zone scattered field 
versus time and the co-polarized RCS versus frequency for the 
dielectric sphere. 

Figures 4-5 show the co-polarized far zone electric field 
versus time and the co-polarized RCS for the magnetic sphere. 

XI. SAMPLE PROBLEM SETUP 

The code as furnished models a 20 cm radius, lossy magnetic 
sphere and computes backscatter far zone scattered fields at 
angles of 0=22.5 and 0=22.5 degrees. The corresponding output 



18 


data files are also provided, along with a code to compute Radar 
Cross Section using these data files. In order to change the 
code to a new problem, many different parameters need to be 
modified. A sample problem setup will now be discussed. 

Suppose that the problem to be studied is RCS backscatter 
versus frequency from a 28 cm by 31 cm perfectly conducting plate 
with a 3 cm dielectric coating with a dielectric constant of 4e n 
using a 0-polarized field. The backscatter angles are 0=30.0 and 
0-60.0 degrees and the frequency range is up to 3 Ghz. 

Since the frequency range is up to 3 Ghz, the cell size must 
be chosen appropriately to resolve the field IN ANY MATERIAL at 
the highest frequency of interest. A general rule is that the 
cell size should be 1/10 of the wavelength at the highest 
frequency of interest. For difficult geometries, 1/20 of a 
wavelength may be necessary. The free space wavelength at 3 GHz 
is A 0 =10 cm and the wavelength in the dielectric coating at 3 GHz 
is 5 cm. The cell size is chosen as l cm, which provides a 
resolution of 5 cells/A in the dielectric coating and 10 cells/A n 
m free space. Numerical studies have shown that choosing the 

cell size < 1/4 of the shortest wavelength in any material is 
the practical lower limit. Thus the cell size of l cm is barely 
adequate. The cell size in the x, y and z directions is set in 
the common file through variables DELX, DELY and DELZ. Next the 
problem space size must be large enough to accomodate the 
scattering object, plus at least a five cell boundary (10 cells 
is more appropriate) on every side of the object to allow for the 
far zone field integration surface. It is advisable for plate 
S f a ^ ering to have plate centered in the x and y directions 
of the problem space in order to reduce the cross-polarized 
backscatter and to position the plate low in the z direction to 
allow strong specular reflections multiple encounters with the 
ORBC. A 10 cell border is chosen, and the problem space size is 
chosen as 49 by 52 by 49 cells in the x, y and z directions 
respectively. As an initial estimate, allow 2048 time steps so 
that energy trapped within the dielectric layer will radiate. 

Thus parameters NX, NY and NZ in COMMONC.FOR would be changed to 
reflect the new problem space size, and parameter NSTOP is 
changed to 2048. if all transients have not been dissipated 
after 2048 time steps, then NSTOP will have to be increased. 
Truncating the time record before all transients have dissipated 
will corrupt frequency domain results. Parameter NZFZ must be 
equal to 1 since we are interested in far zone fields only. 
Parameter MAGNET must be equal to 0 for the dielectric scatterer. 

bU1 ^ t * ie object, the following lines are inserted into the 
BUILD subroutine: 

C 

C BUILD THE DIELECTRIC SLAB FIRST 

C 


ISTART=11 
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JSTART=11 
KSTART=11 
NXWIDE=2 8 
NYWIDE=3 1 
NZWIDE=3 
MTYPE=2 

CALL DCUBE ( ISTART , J START , K START , NXWIDE , NYWIDE , NZWIDE , MTYPE ) 
C 

C BUILD PEC PLATE NEXT 

C 

ISTART=11 

J START= 1 1 

KSTART=11 

NXWIDE=28 

NYWIDE=3 1 

NZWIDE=0 

MTYPE=1 

CALL DCUBE ( I START , JSTART , KSTART , NXWIDE , NYWIDE , NZWIDE , MTYPE) 

The PEC plate is built last on the bottom of the dielectric 
slab to avoid any air gaps between the dielectric material and 
the PEC plate. In the common file, the incidence angles THINC 
and PHINC have to be changed to 30.0 and 60.0 respectively, the 
cell sizes (DELX, DELY, DELZ) are set to 0.01, and the 
polarization is set to ETHINC=1.0 and EPHINC=0.0 for ©-polarized 
fields. Since dielectric material 2 is being used for the 
dielectric coating, the constitutive parameters EPS (2) and 
SIGMA(2) are set to 4e 0 and 0.0 respectively, in subroutine 
SETUP. This completes the code modifications for the sample 
problem. 

XII. NEW PROBLEM CHECKLIST 

This checklist provides a quick reference to determine if 
all parameters have been defined properly for a given scattering 
problem. A reminder when defining quantities within the code: 
use MKS units and specify all angles in degrees. 

COMMONC . FOR : 

1) Is the problem space sized correctly? (NX, NY, NZ) 

2) For near zone fields, is the number of sample points correct? 
(NTEST) 

3) Is parameter NZFZ defined correctly for desired field 
outputs? 

4) Is parameter MAGNET defined correctly for the type of 
scatterer? 

5) Is the number of time steps correct? (NSTOP) 
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6) Are the cell dimensions (DELX, DELY, DELZ) defined correctly? 

7) Are the incidence angles (THINC, PHINC) defined correctly? 

8) Is the polarization of the incident wave defined correctly 
(ETHINC, EPHINC)? 

9) For other than backscatter far zone field computations, are 
the scattering angles set correctly? (THETFZ , PHIFZ) 

SUBROUTINE BUILD: 

1) Is the object completely and correctly specified? 

SUBROUTINE SETUP: 

1) Are the constitutive parameters for each material specified 
correctly? (EPS, XMU, SIGMA, SIGMAC) 

FUNCTIONS SOURCE and DSRCE: 

1) If the Gaussian pulse is not desired, is it commented out and 
the smooth cosine pulse uncommented? 

SUBROUTINE DATSAV : 

1) For near zone fields, are the sampled field types and spatial 
locations correct for each sampling point? (NTYPE, IOBS, JOBS, 
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XIV. FIGURE TITLES 

Fig. 1 Standard three dimensional Yee cell showing placement 
of electric and magnetic fields. 

Fig. 2 Co-polarized far zone scattered field versus time for 
dielectric sphere with 20 cm radius. 

Fig. 3 Co-polarized Radar Cross Section versus frequency for 
dielectric sphere with 20 cm radius. 

Fig. 4 Co-polarized far zone scattered field versus time for 
magnetic sphere with 20 cm radius. 

Fig. 5 Co-polarized Radar Cross Section versus frequency for 
magnetic sphere with 20 cm radius. 
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